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Abstract 

^ ' Our model is a constrained homogeneous random walk in Zi_ . The convergence to stationarity for 

oo : 

\Q . such a random walk can often be checked by constructing a Lyapunov function. The same Lyapunov 

^-j- ■ function can also be used for computing approximately the stationary distribution of this random 

walk, using methods developed by Meyn and Tweedie in p4| . In this paper we show that, for this 
type of random walks, computing the stationary probability exactly is an undecidable problem: no 
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model is also an undecidable problem. We extend these results to a certain type of queueing systems. 
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algorithm can exist to achieve this task. We then prove that computing large deviation rates for this 



The implication of these results is that no useful formulas for computing stationary probabilities and 
large deviations rates can exist in these systems. 

1 Introduction 

The main model considered in this paper is a constrained homogeneous random walk in a <i-dimensional 
nonnegative orthant Zi, where Z+ is the space of d-dimensional vectors with integral nonnegative 
components. Specifically, the transitions with positive probabilities can occur only to neighboring 
states and the transition probabilities depend only on the face that the current state of the random 
walk belongs to, but not on the size of the components of the state. 
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Ever since the appearance of the papers by Malyshev f28|| , [27), |2J| and Menshikov |3l[], random 
walks in Z+ have assumed a prominent role in modelling and analysis of queueing networks of certain 
type, for example Markovian queueing networks. Specifically, the question of positive recurrence or 
stability was analyzed. One of the main techniques used for the stability analysis of these type of 
random walks is a Lyapunov function technique also known as Foster's criteria. A comprehensive 
study of constrained random walks in Z+ was conducted by Fayolle, Malyshev and Menshikov in [15] 



and many additional results appeared after the book was published. Specifically, a very interesting 
connection between constrained random walks and general dynamical systems on compact manifolds was 
established by Malyshev [ 30 1 . Exact conditions for positive recurrence for the case d < 4 were obtained 



by Ignatyuk and Malyshev in [18]. The large deviation principle for special cases and modifications of 
random walks in Z+ was established by Ignatyuk, Malyshev and Scherbakov in |^]. This followed by 
efforts to actually compute the large deviation rates, which turned out to be a very complicated problem. 
See for example Kurkova and Suhov |26|| , where large deviations limits are computed for a random walk 
in Z\ arising from joint-the-shortest queueing system. The analysis uses a fairly complicated complex- 
analytic techniques. The goal of the current paper is to explain the difficulty in obtaining such results 
for general dimensions. 

Analysis of random walks arising from special types of multiclass queueing networks became a subject 
of particularly aggressive research efforts during the previous decade. Many interesting and deep results 
were established which connect stability of such queueing networks with stability of corresponding fluid 
models, obtained by Law of Large Numbers type of rescaling. This research direction was initiated 



in pioneering works by Rybko and Stolyar and Dai ]10|, where it was shown that stability of a 
fluid model implies stability of the underlying queueing system. The converse of this result is not true 
in general, see Dai, Hasenbein and Vande Vate [12], Bramson ||, but is true under some stronger 



conditions, Dai [11], Meyn |32| . Despite these results, to the day no full characterization of stable 
queueing networks is available. Stability was characterized only for special types of queueing networks 
or special scheduling polices. For example feedforward networks are known to be stable for all work- 



conserving policies, Down and Meyn [14], Dai [pLCH . Stability of fluid networks with two processing 



stations operating under arbitrary work-conserving policies is fully characterized in Bertsimas, Gamarnik 



and Tsitsiklis [|l[] by means of a linear programming and in Dai and Vande Vate [13] by direct methods. 
The question of computing stationary distributions comes naturally after the question of stability. 
Several results are available again in the context of multiclass queueing networks. Some of the results 
were obtained using quadratic Lyapunov functions, Bertsimas, Paschalidis and Tsitsiklis M, Kumar 



and Kumar [25], Kumar and Meyn 23], Kumar and Morrison [24], using piece- wise linear Lyapunov 
functions in Bertsimas, Gamarnik and Tsitsiklis Q, and using more direct methods, Bertsimas and 
Nino- Mora ||], 0] All of the results obtain only bounds on the stationary probabilities. Computing 
exactly the stationary probabilities seems beyond the existing techniques. 



2 Our results 

It was established by the author in |[(| that positive recurrence of a constrained homogeneous random 
walk in Z+ is an undecidable property. Meaning, no algorithm can possibly exists which given the 
description of the random walks (given the dimension and the transition matrix) will be able to check 
whether the walk is positive recurrent. This result was also established for queueing systems operating 
under the class of so called generalized 'priority policies. This result explains the difficulty in stability 
analysis by stating that these problems are simply insolvable. It was conjectured in the same paper 
that the stability of multiclass queueing networks operating under the class of much studied priority or 
First-In-First-Out policies, is undecidable as well. The conjecture remains unproven. 

In the current paper we continue the decidability analysis of constrained random walks by asking 
the following question: given a constrained homogeneous random walk, can we compute its stationary 
distribution, provided that the existence of a stationary distribution can be checked, for example, by 
constructing a Lyapunov function? To put this question into a proper computation theoretic framework, 
we ask the following question. Given a constrained homogeneous random walk, which possesses a unique 
stationary distribution it, given a state q € Zi, for example q = 0, and given a rational value r > 0, is 
it true that the stationary probability n{q) of this state satisfies n(q) < r? In this paper we prove that 
this problem is undecidable, even if a Lyapunov function witnessing positive recurrence is available. 
Thus, no algorithm can exist which given a positive recurrent constrained homogeneous random walk 
computes its stationary distribution. Specifically, the stationary distribution cannot be written down in 
any constructive way using some formulas. Contrast this with random walks corresponding to product 
form type networks, for example Jackson networks, for which a very simple formula is available. We 
then prove that computing large deviations limits for the same model is an undecidable problem as 
well. In particular, we show that given a random walk in Z+ with a unique stationary distribution 
(witnessed, for example, by a Lyapunov function), and given a vector v £ $l d , the problem of deciding 
whether lirn^oo \og(-K(vn))/n is finite, is undecidable. We extend these results to queueing systems 
operating under a class of generalized priority policies. Finally, we observe that, nevertheless, estimating 
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stationary distribution of a constrained random walk is a decidable problem, if one is willing to tolerate 
a two-sided error and a Lyapunov function exists. Specifically, given such a random walk with the 
unique stationary distribution 7T, given a Lyapunov function, given a state q £ Z+ and any value e > 0, 
an interval (r, r + e) can be constructed which contains vr(g). This result is an easy consequence of a 



powerful result obtained by Meyn and Tweedie [34|, which obtains exponential bounds on the mixing 



rate of Markov chains, using Lyapunov function methods. We note that such approximation result 
cannot be obtained for large deviations rates since, as we mentioned above, even determining whether 
a given large deviation rate is finite is an undecidable problem. 

The remainder of the paper is organized as follows. In the following section we describe our model 
- constrained homogeneous random walk in Zi and introduce Lyapunov functions. In Section || we 
introduce a counter machine - a modification of a Turing machine which for us is the main tool for 
establishing the undecidability results. In Section || we prove that computing a stationary distribution 
of a positive recurrent random walk in Zi is an undecidable problem. In Section ^ we prove that 
computing large deviations rates for positive recurrent random walks in Zf_ is an undecidable problem 
as well. Extension of these results to queueing systems is established in Section ^. In Section ^ we 
show how stationary distribution can be computed with a two-sided error using a Lyapunov function 



technique and Meyn and Tweedie results from [34]. Conclusions and open problems are discussed in 
Section [| 

3 Constrained homogeneous random walk in Z+. Lyapunov function 
and stationary distribution 

Let Zi denote the space of d-dimensional vectors with nonnegative integer components. Our model is 
a random walk Q(t),t = 0, 1, 2, . . . which has Zi_ as a state-space. For each A C {1, 2, . . . , d} let Z^ 
denote the corresponding face: 

Z A = {( z i , *2, ■ ■ ■ , Zd) G ■ z i > for i £ A, Zi = for i £ A}. 

The transition probabilities are face-homogeneous - they depend entirely on the face the random walk is 
currently on. In addition the transition vectors have at most unit length in max norm. In other words, 
for each A C {1, 2, . . . , d} and each A € {—1, 0, l} d a certain value p(A, A) (the transition probability) 
is defined. These values satisfy 

E P(A,A) = 1 
Ae{-i,o,i} <i 



for each A and p(A, A) = if for some i ^ A, Aj = —1. The latter condition is simply a consistency 
condition which prevents transitions into states with negative components. Given a current state Q(t) G 
i?^ of the random walk, the next state is chosen to be Q(t)+ A with probability p(A, A), if the state Q(i) 
belongs to the face Z\. We will also write p(q, q') instead of p(A, A) if the state q G Z\ and q' — q = A. 
We denote by p^(q,q') the i-step transition probabilities: Prob{<5(t) = q'\Q(0) = q}. The model 
above will be referred to as constrained homogeneous random walk in Zi_. We will say that our walk 
is deterministic if p(A, A) G {0,1} for all A and A. In other words, the transition vector A = A(A) 
deterministically depends on the face. The set of parameters p(A, A) is finite, and, in particular, it 
contains 6 d elements corresponding to 2 d faces Z\ and 3 d transition vectors A per face (with some 
transitions occurring with zero probability). Let ||Q(i)|| denote L\ norm. That is ||Q(t)[| = J2i<dQi(t)- 
For any state q G Zi and subset X C Zi, let T = T(q,X) denote the first hitting time for the set 
X when the initial state of the walk is q, including the possibility T = oo. That is 

T = mm{t : Q(t) G X\Q(0) = q}. (1) 

The following definition is standard in the theory or infinite Markov chains. 

Definition 1 A homogeneous random walk is defined to be positive recurrent or stable if there exist 
some C > such that the random walk visits the set Xq = {z£ Z^_ : Ylt=i z i — C} infinitely often with 
probability one, and E[T(q,Xc)] is finite for all q G Zi_. 

Stability of a constrained homogeneous random walk Q(t) can be checked, for example, by con- 
structing a suitable Lyapunov function. 

Definition 2 A function $ : Z+ — > is defined to be a Lyapunov function with drift —7 < and 
exception set B C if \B\ < 00 and for every state q £ B 

E[$(Q(t + l)\Q(t) =q}- <%) = J2 <%>(^ V') ~ <%) < "7- (2) 

q ez d + 

In other words, the expected value of the Lyapunov function should decrease at each time step, 
whenever the random walk is outside of the exception set. Existence of a Lyapunov function under some 
additional assumptions implies stability. For a comprehensive survey of Lyapunov function methods 
see Meyn and Tweedie |33|. Various forms of Lyapunov functions, specifically piece-wise linear and 
quadratic Lyapunov functions were used to prove stability of random walks corresponding to Markovian 
queueing networks [g], (5j, [21], |23|], [24], pq| . In some simple cases even linear Lyapunov function of 



the form &(q) = J2i=i w iQi with w% > can prove stability of a constrained homogeneous random walk. 
It is easy to see that a linear function <£((/) = w T ■ q is a Lyapunov function if and only if for some 7 > 
and every nonempty set A C {1, 2, . . . , d} the following inequality holds 



E[w T Q(t + 1) - w T Q(t)\Q(t) eZ A ]= {w T A)p(A, A) < -7. (3) 

Ae{-i,o,i} d 

The existence of a linear Lyapunov function is only sufficient but not necessary for stability of the 
constrained random walk. It is also useful sometimes to consider a geometric Lyapunov function, defined 
as follows. 

Definition 3 A function $ 9 : — ► [1, +00) is defined to be a geometric Lyapunov function with drift 
< 7 9 < 1 and exception set B C Zi if \B\ < 00 and for every state q B 



A geometric Lyapunov function is used, for example in Meyn and Tweedie 34], to prove exponentially 
fast mixing of a Markov chain which admits a geometric Lyapunov function. The precise statement of 
this result will be given below in Section || If the condition (||) is met for some w and 7, then a function 
of the form & g (q) = ex.p(5w T ■ q) is a geometric Lyapunov function for a suitable value of 5 > 0. 

Throughout the paper we will assume all the states q communicate with the state 0, that is pC' (q, 0) > 
for some t > 0. As a consequence, the random walk is irreducible. If it is in addition positive recurrent, 
then it possesses a unique stationary distribution ir : Zf — > [0, 1]. Namely X^ez^ 7r (?) = 1 an d f° r an y 
state q 

E t(«W,9)=t(?) (5) 

q'&Z<l 

This stationary distribution is defined completely by the set of transition parameters p(A, A). Com- 
puting the stationary probability distribution for these walks is the main focus of this paper. It was 
established by the author in [16] that checking positive recurrence of a constrained homogeneous random 
walk is an undecidable problem - no algorithm can exist to achieve this task. However, if one is lucky 
to construct a Lyapunov function, for example by checking condition |3| for some nonnegative vector 
w G then the random walk is in fact positive recurrent. One might be tempted to believe that in 
this case the analysis of the random walk is simplified significantly. In Section |5] we show even if a lin- 
ear Lyapunov function exists, computing the stationary probability distribution is still an undecidable 
problem. As in the case of stability analysis, our main tool for establishing this undecidability result is 
a counter machine and the halting problem defined in Section |4|. 



4 Counter Machines, Halting Problem and Undecidability 



A counter machine (see Jl8|]) is a deterministic computing machine which is a simplified version 
of a Turing machine - a general description of an algorithm working on a particular input. In his 
classical work on the Halting Problem, Turing showed that certain decision problems simply cannot 
have a corresponding solving algorithm, and thus are undecidable. For a definition of a Turing machine 



and the Turing Halting Problem see [37]. Ever since many quite natural problems, in mathematics and 
computer science were found to be undecidable. Some of the undecidability results in control theory 
were obtained by reduction from a counter machine, see Blondel et al. ||. For a survey of decidability 
results in control theory area see Blondel and Tsitsiklis [^J. 

A counter machine is described by 2 counters i?2 and a finite collection of states S. Each counter 
contains some nonnegative integer in its register. Depending on the current state s € S and depending 
on whether the content of the registers is positive or zero, the counter machine is updated as follows: 
the current state s is updated to a new state s' € S and one of the counters has its number in the 
register incremented by one, decremented by one or no change in the counters occurs. 

Formally, a counter machine is a pair (S, T). S = {sq, s\, . . . , s m _i} is a finite set of states and T is 
configuration update function r : S x {0, l} 2 — > S x {—2, —1, 0, 1, 2}. A configuration of a counter ma- 
chine is an arbitrary triplet (s, z\, Z2) £ S x Z\. A configuration (s, Z\, Z2) is updated to a configuration 
(s',z'i,Z2) as follows. First a binary vector b = (61,62) is computed were bi = 1 if Zi > and hi = if 
Zi = 0, i = 1, 2. If r(s, b) = (s', 1), then the current state is changed from s to s', the content of the first 
counter is incremented by one and the second counter does not change: z'± = Z\ + 1, ^ = %%■ We will 
also write T : (s, z%, Z2) —* (s' , Z\ + 1, Z2) and r : s — * s', T : z\ —* z\ + 1, T : Z2 — ► Z2- If T(s, b) = (s', —1), 
then the current state becomes s' , z[ = z\ — 1, z' 2 = Z2- Similarly, if T(s, b) = (s', 2) or T(s, b) = (s', —2), 
the new configuration becomes (s', z%, Z2 + 1) or (s', zi, Z2 — 1), respectively. If T(s, b) = (s', 0) then the 
state is updated to s', but the contents of the counters do not change. This definition can be extended 
to the one which incorporates more than two counters, but, in most cases, such an extension is not 
necessary for our purposes. 

Given an initial configuration (s°, Zj 1 , z!>) the counter machine uniquely determines subsequent con- 
figurations (s , z\, zl), {s 2 , z 2 , z 2 ), . . . , (s*, z[, Z2), ... . We fix a certain configuration (s*, zl, z^) and call 
it a halting configuration. If this configuration is reached then the process halts and no additional 
updates are executed. The following theorem establishes the undecidability of the halting property. 
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Theorem 1 Given a counter machine (S,T), initial configuration (s , Zi,z%) and the halting configu- 
ration (s*,^,2|), the problem of determining whether the halting configuration is reached in finite time 
is undecidable. It remains undecidable even if the initial and the halting configurations are the same 
with both counters equal to zero: s° = s*, = z\ = z* = z\ = 0. 



The first part of this theorem is a classical result and can be founded in [17]. The restricted case of 
s° = s*, z® = z*,i = 1,2 can be proven similarly by extending the set of states and the set of transition 
rules. It is the restricted case of the theorem which will be used in the current paper. 



5 Computing the stationary probability distribution. The undecid- 
ability result 

Theorem |l| was used in |ll| to prove that the stability of a constrained random walk in Zf_ is undecidable. 
Naturally, the problem of stability comes before the problem of computing the stationary distribution of 
a stable random walk. As we mentioned in Section ^, stability can be checked sometimes by constructing 
a Lyapunov function. In this section we prove our main result: even if such a Lyapunov function, 
witnessing stability, is available and is provided as a part of the data parameters, computing stationary 
distribution is an undecidable problem. 

We now give an informal outline of the proof. The proof uses a reduction from a halting problem for 
a counter machine. We embed a counter machine with initial and halting configuration (s*,0, 0) into 
a deterministic walk in Z+. The state space and the transition rules of this walk are then extended in 
some way that incorporates an independent Bernoulli process with some fixed parameter p. We then 
show that 

• If the original counter machine never returns to the initial configuration (s*,0, 0), then the con- 
structed random walk, when started from the origin, returns into the origin in 2t + 2 steps with 
probability (1 — p)p for i = 0,1,2,... .In particular, the expected return time to the origin is 
2/(1 -p). 

• If the original counter machine returns to the initial configuration in T steps, then the modified 
random walk returns into the origin in 2t + 2 steps with probability (1 — p)p l for t < T — 1 and in 
2T + 2 steps with the remaining probability 1 — J2t<T-i (1 ~p)p l = P T ■ I n particular, the expected 
return time to the origin is (2 — 2p T+l )/(l —p). 
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The stationary probability distribution of any state is exactly the expected return time to this state. 
Therefore, the stationary probability of the origin is (1 — p)/2 if the counter machine halts and is 
strictly greater, if the counter machine does not halt. Since the value p is our control parameter, and 
since checking whether the counter machine halts is an undecidable problem, then computing stationary 
probability is undecidable as well. We now state and prove rigorously this result. As before, let tt denote 
the unique stationary distribution of an irreducible positive recurrent random walk. Let also denote 
the origin of the nonnegative lattice Zi. . 

Theorem 2 Given an irreducible constrained random walk with transition probabilities p(A, A), given 
a linear vector w G satisfying ^) and given a rational value < r < 1, the problem of checking 
whether vr(0) < r is undecidable - no algorithm exists which achieves this task. 

Remarks : 

1. The stationary distribution can in principle take non-rational values. In order to put the problem 
into a framework suitable for algorithmic analysis we modified the question into the one of checking 
whether ir(-) < r for rational values r. This is a standard method in the theory of Turing decidable 
numbers, see g. 

2. A simple example where computing the stationary probability distribution is a decidable problem 
is Jackson networks, For such a network with d stations the stationary probability of the state 
m = (mi, . . . , md) is given by TJ^ =1 (1 — Pj)p^ 3 , where pj is the traffic intensity in station j. Specifically, 
the stationary probability of the state is rX/ = i(l ~~ Pj)- Given any rational value < r < 1, it is a 
trivial computation to check whether this product is at least r. 

Proof of Theorem ^: we start with a construction used in [jl6|. Namely, we embed a given counter 
machine with states Sq, s\, . . . , s m _i into a deterministic walk in as follows. Without the loss of 

generality, assume that s* = so- Let configuration (sj, z\, Z2), 1 < i < m — 1 correspond to the state 
q = (ej, ZijZz) G where is unit vector with 1 in i-th coordinate and zero everywhere else. Also, 

let configurations (so, z±, Z2) correspond to (0, z\, Z2), with zeros in first m — 1 coordinates. Specifically, 
the initial and halting configuration (so,0, 0) corresponds to the origin 0. We now describe the set of 
transition vectors A = A (A). We describe it first for subsets A C {1,2, . . . ,m + 1} which correspond 
to an encoding of some configuration of a counter machine. Specifically, A n {1,2, . . . ,m — 1} = 
(corresponding to configurations with state sq) or An {1,2, ... ,m — 1} = {i} for some 1 < i < 
m — 1, corresponding to configurations with state S{. Fix any configuration (si,zi,Zz). Suppose the 
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corresponding update rule is r((sj, z±, zz)) = (sj,+l) for some 1 < j < m — 1. That is, the state is 
changed into Sj, the first counter is incremented by 1 and second counter remains unchanged. We make 
the corresponding transition vector to be A = A (A), where the i-th coordinate of A is —1, the j-th 
coordinate is +1, the m-th coordinate is +1 and all the other coordinates are zeros. It is easy to see 
that if at time t, the state Q(t) corresponds to some configuration (s^ Zx,z%), that is Q(t) = (ej, £i, #2), 
then Q(t + 1) = Q(t) + A corresponds to the configuration (sj, z\ + 1, z-z) obtained by applying rule T. 
We construct transition vectors similarly for other cases of configuration updates. In particular, if the 
state Si is changed to state sq, then the corresponding A has —1 in the i-th coordinate and zeros in all 
the coordinates 1 < j < m — 1, j ^ i. As we will see later, if Q(t) corresponds to some configuration of 
a counter machine at time t, then it does so for all the later time t' > t. Now if Q(t) belongs to some 
face 2\ which does not correspond to some configuration, then we simply set A(A) = — e» where i is the 
smallest coordinate which belongs to A. Then at some later time t' > t the state Q{t') will correspond 
to some configuration. 



Construction above is exactly the one used in [16| to analyze stability. We now modify the con- 



struction by adding two additional coordinates. Our new state at time t is thus denoted by Q(t) = 
{Q(t),q\(t),q2(t)) € iT^" 1 " 3 . Also a parameter < p < 1 is fixed. The transition rules are modified as 
follows. 

1. When q2{t) = 1, the first part Q(t) of the state is updated exactly as above. Also, if > 0, in 
other words, Q(t) does not represent the halting configuration (so, 0, 0), then the value of <72(i) stays 
1 with probability p and switches to with probability 1 — p. If, on the other hand ||Q(i)|| = 
then we set 52 (i + 1) = with probability 1. Finally, the value of qi(t) £ {—2,-1,0,1,2} is 
selected in such a way that \\(Q(t + l),qi(t + 1))[| = ||(Q(i), + 1, where = 
Si^i 1 Qi(t) + Qi(t)- It is eas Y to see that such a value of q± (t) always exists. For example if Q(t) 
encodes (sj, z±, z%),i 7^ and the configuration is changed into (sj, Z\,Z2 — 1), j 7^ 0, then we put 
qi (t) = 2. 

Remark : We stipulated before that the transition vectors A must belong to {— 1,0, l} m+4 for 
our constrained random walk, whereas above the value of q\{t) can change by —2 and 2. It is 
easy to satisfy this constraint for q\{t) by splitting it into two coordinates and making 

qi (t) = q[(t) = 1 in case q±(t) was assigned 2 before, and q±(t) = q[(t) = —1 in case qi(t) was 
assigned —2. We keep only one qi(t) for simplicity, allowing it to take values —2,2. 

2. When q^if) = 0, we set A^ = —1, Aj = 0, i 7^ k, 1 < i < m + 3, where k is the smallest coordinate 
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such that Qk(t) > 0. In particular, q2(t) stays equal to 0. If Q(t) = qi(t) = 0, (in particular 
Q(t) encodes the initial-terminal configuration (so, 0,0)) then Q(t) and qi(t) are updated as in 
the case q% (t) = 1 above. Also q%{£) in this case is switched to 1 with probability p and stays 
with probability 1 — p. 

Note, that the only stochastic part in our random walk is the last component q^if). 

Proposition 1 The constructed random walk Q(t) is irreducible and positive recurrent with the unique 
stationary distribution it. Moreover, 

1. If the counter machine with the initial configuration (sq,0, 0) does not halt, then the random walk 
Q(t) with the initial state 0(0) = returns to the origin in2t + 2 steps with probability (1 — p)p t , 
for t = 0, 1, 2, . . . . As a result, the expected recurrence time of the state is l/ir(Q) = 2/(1 — p). 

2. If the counter machine with the initial configuration (so, 0, 0) halts inT > 1 steps, then the random 
walk Q(t) with the initial state 0(0) = returns to the origin in 2 + 2t steps with probability 
(1 — p)p t for t < T, and in 2 + 2T steps with the remaining probability p T . As a result, the 
expected recurrence time of the state is l/vr(0) = (2 — 2p T+1 )/{l — p). 

3. For any C > 2/(1 -p) the function J2Zf Qi + + Cq 2 (t) is a linear Lyapunov function with 
drift —7 = — 1 and an exception set B = {0}. 

We first show that the proposition above implies the theorem. Suppose, we had an algorithm A 
which given an irreducible constrained random walk Q(t), with a linear Lyapunov function w T Q(t) 
and given a rational value < r < 1 could determine whether the unique stationary distribution ir 
satisfies vr(0) < r. We take a counter machine and construct a random walk Q(t) as described above. 
Proposition implies that this walk is a valid input for the algorithm A. We use A to determine 
whether vr(0) < r = (1 —p)/2. From Proposition |l], this is the case if and only if the underlying counter 
machine does not halt. In this fashion, we obtain an algorithm for checking halting property for counter 
machines. This is a contradiction to Theorem |l|. □ 

Proof of Proposition Suppose the underlying counter machine does not halt. Let us trace the 
dynamics of our random walk Q(t) starting from Q(0) = 0. Initially, by applying rule 2, it moves into 
some state (Q(1),0, 1) with probability p or state (Q(1),0,0) with probability 1 — p. An independent 
Bernoulli process for q^ify with parameter p is continued in the first case. Suppose this process succeeds 
exactly t > times (including the transition from initial state 0), which occurs with probability (1— p)p t . 
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Then, applying rule 1, at times t and t + 1 we have states (Q(t), qi(t), 1), (Q(t + l),qi(t + 1),0) with 
\\Q(t) + 9i (01 1 = \\Q(t + 1) + Qi(t + 1)|| = t + 1. At this moment rule 2 becomes applicable. Since at 
each step the norm ||Q(i)|| decreases exactly by one, the origin is reached at time (t + 1) + (t + 1)- We 
conclude that the return time is 2 + 2t with probability (1 — p)p t ,t = 0, 1, 2, . . . . The expected return 
time is then 2/(1 — p) and the stationary probability of the state is (1 — p)/2. 

Suppose, now, the underlying counter machine reaches the terminal state (so,0, 0) in exactly T > 1 
steps. Suppose also the Bernoulli process for q2(t) succeeds exactly t > times for t < T. Then, 
exactly as above, the origin is reached in 2 + 2t steps and this occurs with probability (1 — p)p . If, 
however, by the time T the Bernoulli process does not fail, which occurs with probability p T , then the 
state Q(T) = (Q(T),qi(t),l) = (0 m +i, <7i(i), 1) is reached at time T, where m+ i denotes a m + 1- 
dimensional zero vector. By the choice of qi(t) in rule 1, qi(T) = T. At time T + 1, by rule 1, we have 
a state Q(T + 1) = (0 m+ 2, T + 1, 0) and rule 2 applies. At time (T + 1) + (T + 1) the origin is reached. 
We conclude that the random walk returns to the origin in 2 + 2T steps. Combining, the expected 
return time to the origin is then (2 — 2p T+1 )/(l — p), if the counter machine halts in T steps, and the 
stationary probability of the state is (1 — p)/(2 — 2p T+1 ). 

To complete the proof of the proposition, we analyze the expected change of the function Y^ 1 i=\ L Qi + 
qi(t) + Cq2(t) . When q2(t) = and Q(t) ^ 0, the sum decreases deterministically by 1. When q2(t) = 1, 
the value of X^l* Qi + 9l(*) increases deterministically by 1, and the value of Cq^it) stays the same 
with probability p or decreases by C with probability 1 — p. Therefore, the expected change of the sum 
is 1 — C(l — p). When C > 2/(1 — p), the expected change is at most —1. □ 

An important implication of Theorem [2] is that it is impossible to express the stationary distribution 
of a positive recurrent random walk Q(t) as a function of the parameters p(A, A) via some computable 
function /(•)• For example, the stationary distribution cannot be expressed as roots of some polynomial 
equations with rational coefficient, as inequalities x < r can be checked for any root x of such a 
polynomial and any rational value r. This is a startling contrast to a simple expression Y\(l — Pj)Pj j 
corresponding to a stationary distribution of a Jackson network. 

6 Large Deviation Rates. The undecidability result 

In this section we discuss the question of computing large deviation rates for our model. Specifically, we 
focus on computing large deviation rates for the stationary distribution tt of our random walk Q(t) in 
Zf. Let 3? and !R + denote the set of real values and the set of nonnegative real values, respectively. For 
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any x G 5ft let \x\ denote largest integer not bigger than x, and for any x E 5R d let \x\ = ( \x\\ , • • • , \xd\ )• 
We say that a function L : — > U {00} is a large deviation rate function for a given irreducible 
positive recurrent random walk Q{t) in .2^ if for any vector v € the stationary distribution it 
satisfies 

lim lo g(vr(M)) =L 

In other words, the stationary probability of being in state [vn\ is asymptotically exp(— L{v)n) 
for large n. For results on large deviations for specific types of constrained random walks in Zi see 
p0|]. There are numerous works on large deviation in the context of queueing systems, see Shwartz 



and Weiss for a survey. Specifically, Kurkova and Suhov [26| study large deviation rates for a two 
dimensional random walk corresponding to join-the-shortest queue. The analysis is quite intricate and 
uses complex-analytic techniques developed by Malyshev p8|] , [p?|], p9| back in 70's. To the best of our 
knowledge, the existence of the large deviations limits (^) is not fully proved for general constrained 
homogeneous random walks Q(t) in Z+. One can instead consider limits 

t r \ v ■ t logML 1 "*])) r < \-v log(vr(LunJ)) 

Liv) = hmmf , L + (v) = hmsup . (7) 

n->oo n n-,00 n 

The goal of the present section is to prove that computing the large deviation rate function L(v) is 
an undecidable problem, even if the walk is known to be a priori positive recurrent via, for example, 
existence of a linear Lyapunov function, and even if the large deviation limit function L(v) is known to 
exist. The following is the main result of this section. 

Theorem 3 Given an irreducible constrained random walk, given a linear vector w G satisfying 
^j, given a rational value < r < 1 and a vector v € Z^, the problems of determining whether 
L-(v) < r,L + (v) < r are undecidable. 

Remark : As we will see below, the large deviations limit function L(v) = L^(v) = L + (v) exists 
for the subclass of random walks we consider. As before, the reason for including a linear Lyapunov 
function into the condition of the theorem is to provide a simple way of insuring that the walk is positive 
recurrent. 

Proof : The proof is again based on reduction from a halting problem for a counter machine. Given 
a counter machine with m states consider the extended m+3-dimensional random walk Q(t) constructed 
in the proof of Theorem [2|. We extend it even further by adding an additional coordinate (Q(t), qs(t)). 

13 



Recall, that the rules for updating q\{t) were such that ||Q(i)|| + qi(t) = t as long as q2 (f) = 1 for 
1 < tf < t. Construct the rules for updating qs(t) as follows. If 32 (i) = 1, then qa(t + 1) = 33 (i) + 1. 
Also, if Q(t) = qi(t) = q2(t) = then again q$(t + 1) = 33 (i) + 1. In other words, as long as the random 
walk starts from the origin and as long as q2{t) remains equal to 1, qz(t) = HQWII + Qi(t) = t- Once 
q2(t) becomes zero, the value of q% (t) stays the same as long as Q(t) 7^ and decreases by one when 
Q(t) becomes zero, and continues decreasing until it itself becomes zero. 

Let v = (0, ... ,0, 1) be an m + 4 dimensional vector with the last coordinate equal to unity and 
all other coordinates equal to zero. We now analyze the large deviation rate L{v) for this vector with 
respect to the unique stationary distribution tt. Specifically, we show that the value of L{v) depends on 
whether the counter machine halts. Indeed, if the counter machine halts in T steps, then the value of 
||Q(*)I| + Qi(t) + Q2(t) + Q3(t) is bounded by 2T + 1 and, as a result, the stationary probability of the 
state nv, for large n becomes zero. That is L(v) = +00. Now we show that if the counter machine does 
not halt, then L(v) = log p. We compute ir(nv) by computing the expected return time to state nv, 
when the random walk is in this state at time 0. Thus, we have qs(0) = n,Q(0) = qi(0) = 92(0) = 0. 
By the update rules of q-j, it decreases by one at each time step and at time t = n it becomes zero. All 
the other components remain equal to zero. Beginning from this time t = n, the random walks keeps 
returning to the origin after some random time intervals. It is easy to see that the probability that the 
state nv is visited in between any given two visits to the origin is exactly p n - the probability that the 
Bernoulli process survives at least n steps. Let i?i,i?2, • • • denote the random time intervals between 
successive visits to the origin. For a fixed m > 1 the probability that R m is the first interval during 
which vn is visited is p n [l — p n ] m ~ l , m = 1,2,..., and the expected number of intervals R m before state 
nv is visited for the first time is (p™) -1 . Let I ri denote the indicator function for the event "state vn 
is visited between visits to the origin". In particular, Prob{/ n } = p n . We now compute E[R m \I n ] and 
-E[i? m |-T n ]. Note, that the time n which takes to get from nv to the origin plus the expected time it 
takes to get from the origin to nv, conditioned on I n is exactly E[R m \I n ]. We then obtain that the 
expected recurrence time of the state nv is {p n )~ l E[R m \I n ] + E[R m \I n ]. 

To compute -E[ii m |-], recall from Proposition [jjthat 

E[R m ] = £[iU/„]Prob{/ n } + E[R m \T n ]Piob{T n } = 2/(1 - p). (8) 
If the Bernoulli process survives t > n steps then R m = 3t + 1 steps. Then 
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E[R m I n ] E f >n(l-Py(3t + 1) (4-p)p" 

= p^iu = — ^ = = A ~ P (9) 

Then we obtain from (0) 

2/(l-p)-(4- P )p n 



1 _ p n 

We conclude that the expected return time to the state nv is 

1 - OT -. »/d-rt-(4-^ 4 _ p)< 



ir(nv) 1 — p n 

and 

log(vr(un)) 
hm = logp 

TO — >oo 72 

as we claimed. We see that the value of L{v) depends on whether the underlying counter machine halts 
or not. Specifically, by taking any rational value r > logp, we conclude that the problem of checking 
whether L(v) < r is undecidable, by appealing again to Theorem [l]. □ 
Remarks : 

1. Note that we cannot determine the value of L(v) even approximately, as we cannot distinguish 
between the cases L(v) < +oo and L(v) = +oo. Contrast this with the results of Section |8[ 

2. We would not need extra coordinate 53 (i) if we were interested in large deviation rate 
linin^oo 7r(||Q(i)ll)/ n °f the stationary distribution of the norm of the state. The analysis would be 
identical to the one above. 

7 Application to queueing systems 

The results of the previous sections have implications to a certain type of queueing systems. A queueing 
system consisting of a single station processor and operating under a certain class of generalized priority 



policies was introduced in 16]. It was shown that, similar to constrained random walks, determining 
stability for these queueing systems is an undecidable problem. In this section we consider the same class 
of system and show that computing stationary probabilities and large deviation rates are undecidable 
problems as well. 

We start with the description of the system. Consider a single station queueing system Q consisting 
of a single server and I types of parts arriving externally. The parts corresponding to type i = 1,2, ... , I 
visit the station Jj times. On each visit each part must receive a service before proceeding to the next 
visit. Only one part among all the types can receive service at a time. While waiting for service for the 
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j-th time, the type i part is stored in buffer Bij. We denote by n the total number of buffers n = J2i=i Ji- 
The service time for each part in each visit is assumed to be equal to unity. Each part can arrive into 
the system only in times which are multiples of some fixed integer value M. Specifically, certain values 
< Pi < 1 are fixed for each type i. For each type i and each m = 0, 1, 2, . . . , exactly one part arrives 
at time mM with probability pi and no part arrives with probability 1—pi, independently for all m and 
all other types. In particular, interarrival times are geometrically distributed with expected interarrival 
time equal to 1/Aj = M/pi, where, correspondingly, Aj is the arrival rate for type i. 

A scheduling policy u is defined to be a generalized priority policy if it operates in the following 
manner. A function u : {0, l} n — > {0, 1, 2, ... , n} is fixed. At each time t = 0, 1, 2, . . . the scheduler looks 
at the system and computes the binary vector b = (pi, 62, ■ ■ ■ , b n ) G {0, l} n , where bi = 1 if there are 
parts in the i-th. buffer and bi = 0, otherwise. Then the value k = u(b),0 < k < n is computed. If A; > 
then the station processes a part in the A>th buffer. If k = the server idles. The map u is assumed 
to satisfy the natural consistency condition: u(b) = k > only when 6& = 1. That is, processing can 
be done in buffer k only when there are jobs in buffer k. Note that the generalized priority scheduling 
policy is defined in finitely many terms and is completely state dependent - the scheduling decision at 
time t does not depend on the state of the queueing system at times t' < t. A usual priority policy 
corresponds to the case when there is some permutation 9 of the buffers {1,2, ... ,n} and u(b) = k 
if and only if = 1 and bi = for all i such that 9(i) < 9(k). In words, priority scheduling policy 
processes parts in buffers with lowest value (highest priority) 9, which still has parts. Once we specify 
the queueing system Q and some generalized priority policy u we have specified some discrete time 
discrete space stochastic process. This process considered in times t = mM, m = 0, 1, 2, . . . is in fact a 
Markov chain. 

Given a generalized priority policy u, a pair (Q,u) is defined to be stable if there exists a finite 
number C > such that the total number of parts in the queueing system Q at time t does not exceed 
C for infinitely many t with probability 1. In other words, the underlying Markov chain is positive 
recurrent. In this case there exists at least one stationary probability distribution. It is known that the 
necessary condition for stability is the following load condition 

/ Ji 

^EE A «< L ( 10 ) 

i=i j=i 

This condition is also sufficient for stability if the policy is work conserving, which does not apply 
here, since we allow idling u(b) = 0. We assume that the load condition above holds. We define a 
Lyapunov function and large deviations rates L(v) for this queueing system in the same way we did for 
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constrained homogeneous random walks in Section ^. As for constrained random walks, we now show 
that computing stationary probability distributions and computing large deviations rates for queueing 
systems operating under generalized priority policies is not possible. As for constrained random walks, 
we show that these problems are impossible to solve even if the underlying Markov chain is known to be 
irreducible and a linear Lyapunov function is available. Let tt denote the unique stationary distribution 
of a given irreducible positive recurrent queueing system (Q,u). Let also denote the state of the 
system with all buffers empty. 

Theorem 4 Given a queueing system Q operating under some generalized priority policy u, given a 
linear Lyapunov function <I> and given a rational value < r < 1, the problem of determining whether 
7r(0) < r is undecidable. Likewise, given a vector v, the problem of determining whether L(v) < r is 
undecidable. 



Proof : A reduction from a counter machine to a queueing system operating under some generalized 
priority policy was constructed in This reduction had the following features. Given a counter 

machine with m states, the corresponding queueing system had 24 buffers and / = 3m + 7 arrival 
streams. There is a one-to-one correspondence between the configurations of the counter machine and 
states of the queueing system. In particular, if a counter machine has configuration (si,z\,Z2) at time 
t, then the queueing system at time (3m + 26)i is in a state which corresponds to this configuration 
in some well-defined way. We omit the details of this reduction and instead refer the reader to |I~6| ]. 
We now modify the reduction to incorporate the extended random walk Q(t) = (Q(t),qi(t),q2(t)) that 
was constructed in Section ||. Recall, that the part Q(t) of this walk represented exactly m states and 
the two counters of the underlying counter machine. We add two additional streams of arrivals which 
correspond to coordinates qi(t) and q2{t). We also construct additional buffers for q\ and q<i exactly 
in the way we did in for counters Z\,Z2- The interarrival times for all the arrival streams, except 
for the stream corresponding to q2, are deterministic and equal to some integer M which is selected to 
be bigger than the number of buffers. For the stream corresponding to q2, at most one part arrives at 
times Mt, t = 0, 1, 2, . . . independently for all t, and the probability that a part does arrive at time t is 
equal to p, where p is the parameter selected in construction of the random walk Q(t). Thus, Pi = p for 
the arrival streams corresponding to q2 and Pi = 1 for all the other arrival streams. Finally, we modify 
the rules of the generalized priority policy to incorporate the rules by which the values of q\ (t) , q2 (t) are 
updated. This can be done in a way similar to the rules corresponding to z±, Z2 in [16|. We thus obtain 
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a system which mimics the dynamics of Q(t) at times Mt, t = 0, 1, 2, . . . . A linear Lyapunov function 
can be constructed again, provided that the parameter p is sufficiently small. Arguing as in the proof 
of Theorem |2], we show that the problem of checking whether 7r(0) < r is undecidable. Similarly, we 
show that the problem of checking whether L(v) < r is undecidable, where v is the unit vector with 
one in the coordinate corresponding to q2{t) and zero in all the other coordinates. For the latter case of 
computing large deviations rates, we add an additional arrival stream and buffers to represent the part 
Q3- □ 



8 Computing stationary probabilities approximately using a Lya- 
punov function 

In this section we show that, despite the results of Section [5], computing the stationary probability 
is possible, if we are willing to tolerate some two-sided error and a computable geometric Lyapunov 
function & g exists. Our result is a simple consequence of the following result established by Meyn and 



Twccdic [34], which shows that infinite state Markov chain mixes exponentially fast when a geometric 



Lyapunov function can be constructed. The following is Theorem 2.3 proven in [34]. 

Theorem 5 Given an irreducible Markov chain Q(t), suppose 3> 9 is a geometric Lyapunov function 
with a geometric drift 7 5 < 1 and the exception set B. Suppose also that tt is the unique stationary 
distribution. Then, there exist constants R > 0,0 < p < 1 such that for any state x E X and any 
function (f) : X — > 5? satisfying 4>{x) < Q{x),M x € X, the following bound holds 

| Y, <KV) (Prob{Q(t) = y\Q(0) = x} - vr(y)) | < %(x)Rp t . (11) 
The constants R,p are computable functions which depend on j g , max^gg $ g (x) and 

(x r ) 

uf = max { 8 :p(x,x') > 0}, p£ in = mm p(x,y). (12) 

Exact formulas for computing R, p are provided in |?4j] . They are quite lengthy and we do not repeat 
them here. These formulas give meaningful bounds only in case < 7 g < 1; Ug < oo;p^ in > 0. 

Given a fixed state xq £ consider the function 4>{xq) = l/$ g (xo),4>(x) = 0,x ^ xq. This 
function satisfies the conditions of the theorem and one obtains a computable bound on the difference 
|Prob{(5(i) = xq\Q(Q) = x} — tt(xq)\, which decreases exponentially fast with t. This bound can be used 
for computing stationary probability distribution tt. 
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Theorem 6 Given a constrained random walk Q(t) in Zi_, given a state xq € and an arbitrary value 
e > 0, under the conditions of Theorem there exists an computable value x which satisfies ir(x) € 
[x — e,x + e]. In other words, the stationary probability of the state xq can be computed approximately 
with an arbitrary degree of accuracy. 

Proof : The proof is a simple consequence of Theorem [5|. We fix an arbitrary initial state Q(0), 
say Q(0) = 0. Compute the values R, p, l/Q g (xo). Select t large enough, so that ^ g (Q(0))Rp t < e. 
Compute the transient probability Q(t) = xq conditioned on Q(0) = 0. This can be done by direct 
calculation since t is finite and from any state there are only finitely many neighboring states that can 
be entered with positive probability. The value Prob{Q(t) = xo\Q(0) = 0} can be taken as x, using 
inequality (|i"l|) and by the choice of t. □ 

As we mentioned above, a similar result cannot be established for large deviations rates L(v), since 
the value of L{v ) changes between logp and +oo depending on whether the underlying counter machine 
halts or not. Therefore, computing the value of L{v ) even approximately still is an undecidable problem. 

9 Conclusions 

We considered in this paper the problems of computing stationary probability distributions and large 
deviations rates for constrained homogeneous random walks in Z^_. Both problems were shown to be 
undecidable - no algorithmic procedure for solving these problems can exist. An implication of these 
results is that no useful formulas for computing these quantities, for example along the lines of formulas 
for product form networks, can exist. For the problems of computing stationary probabilities, we showed 
that an approximate computation is possible with arbitrary degree of accuracy if a suitable geometric 
Lyapunov function can be constructed. Yet the problem of computing large deviation rates remains to 
be undecidable even in approximation sense as even checking whether a large deviation rate along a 
given vector is finite or not, is an undecidable problems. 

We conjecture that these problems remain to be undecidable in more restrictive and interesting class 
of Markov chains corresponding to multiclass queueing networks operating under more conventional 
scheduling policies like First-In-First-Out or priority polices. 
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